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Abstract 

To definite and compute differential invariants, like curvatures, for 
triangular meshes (or polyhedral surfaces) is a key problem in CAGD 
and the computer vision. The Gaussian curvature and the mean cur- 
vature are determined by the differential of the Gauss map of the 
underlying surface. The Gauss map assigns to each point in the sur- 
face the unit normal vector of the tangent plane to the surface at this 
point. We follow the ideas developed in Chen and Wu ,2] (2004) and 
Wu, Chen and Chim(2005) to describe a new and simple approach 
to estimate the differential of the Gauss map and curvatures from the 
viewpoint of the gradient and the centroid weights. This will give us 
a much better estimation of curvatures than Taubin's algorithm |10j 
(1995). 

1 Introduction 



The tensors of curvatures on a regular surface S in the 3D Euclidean space are 
important differential invariants in the theory of surfaces and its applications. 
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Gaussian curvatures and mean curvatures are determined by the differential 
of the Gauss map on the surface S. Given a basis on the tangent plane of 
a point p on the surface E, the differential dN of the Gauss map N can be 
realized as a 2 x 2 symmetric matrix A. In fact, the Gaussian curvature and 
the mean curvature at the point p can be computed from the determinant 
and the trace of the matrix A. Since 1990, many methods, like Chen and 
Schmitt [I] (1992) and Taubin [TUj (1995), to estimate these curvatures are 
proposed. However, most of them are devoted to the investigation of the 
principal curvatures, but not directly from the differential of the Gauss map. 
Usually, the accurate estimation of curvatures at vertices on a triangular 
mesh plays as the first step for many applications such as simplification, 
smoothing, subdivision, visualization and registration, etc. 

Chen and Schmitt pOj (1992) and Taubin ^^(1995) employed the circular 
arcs to approximate the normal curvatures. Their methods need to estimate 
the Euler formula and may cause large errors. In this note, we follow the 
methods developed by Chen and Wu 013(2004, 2005) and Wu, Chen and 
Chi jTTj (2005) to estimate the differential dN of the Gauss map directly and 
then obtain accurate estimations of the Gaussian curvature and the mean 
curvature. This method follows the line of the differential of the Gauss 
map more directly and thus provides us a conceptually simple algorithm to 
estimate curvatures on a triangular mesh. Since our method is more natural, 
the estimation turns out to be more accurate. Indeed, it performs much 
better than many other proposed methods. In section two we recall the basic 
theory about the differential of Gauss map on a regular surface. In section 
three, we briefly review the methods of Chen and Schmitt, and Taubin to 
estimate curvatures on triangular meshes. In section four, we present our 
method for estimating the differential of the Gauss map and curvatures on a 
triangular mesh. In section five, we compare the results of our method with 
Taubin's method. 

2 The Gauss map and curvatures on regular 
surfaces 

Consider a parameterization x : U — * S of a regular surface S at a point p, 
where U is an open subset of the 2D Euclidean space M^. We can choose, at 
each point q of x{U), an unit normal vector N{q). The map : x{U) S"^ is 
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the local Gauss map from an open subset of the regular surface S to the unit 
sphere S"^ in the 3D Euclidean space M^. The Gauss map is differentiable 
and its differential dN of at p in S is a linear map from the tangent space 
into itself. 

Given an orthogonal basis {ei, 62} for the tangent space TSp, we can find 
a 2 X 2 matrix A to represent the Gauss map as follows: set aij to be the 
inner product of the vector Cj with dNp{ej). That is, 

ttij =< ei,dNp{ej) > . (1) 

Thus the matrix A = {aij) represents the linear map dNp. Given a vector 
V = aiCi + Q!2e2 in TSp, dNp{v) is equal to w = (3iei + /32e2 where j3i and (32 
satisfy 



(2) 



' I3i\ ^ f ail (ii2\ f ai 

Under this representation of the linear map dNp, the Gaussian curvature 
K and the mean curvature H can be computed hy K = det(A) and H = 
__trace(A)^ Namefy 

K = aiia22 — «i2a2i (3) 
^^_an + a22 

Indeed, since the linear map dNp is self-adjoint, the matrix A should be 
symmetric and diagonalizable. Its eigenvalues are —ni and —K2. The values 
Ki and K2 are the principal curvatures and their associated eigenvectors vi 
and V2 are called the principal directions. In terms of the principal curvatures 
Ki and K,2, we have the Gaussian curvature K = K1K2 and the mean curvature 

TT _ K1+K2 
11 2 • 

Consider a regular curve c(s) with arc length in the regular surface S and 
c(0) = p. The number /c„ =< c"(0), A^(p) > is called the normal curvature 
of c(s) at p. Meusnier's Theorem |S](1976) implies that all curves at p with 
the same tangent vector will have the same normal curvature. This allows 
us to speak of the normal curvature along a unit tangent vector at p. The 
maximum and the minimum normal curvatures are nothing but the principal 
curvatures Ki and K2- Moreover, we have, for any unit vector v in TSp, 

V = vi cos 6 + V2 sin 6 (5) 
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where Vi and V2 are the principal directions and 9 is the angle between 
the vectors v and Vi. The normal curvature /c„ along the unit vector v is 
given by 

kniv) = f^i cos^ 9 + K2 sin^ 9 (6) 
This is known as the Euler formula. 

3 The normal curvature approach: Chen and 
Schmitt's method and Taubin's method 

In this section we shall review the methods given by Chen and Schmitt 
[lj(1992) and Taubin JU](1995). First we introduce some notations. Con- 
sider a triangular mesh S = {V,F), where V = {vi\l < i < n^} is the list 
of vertices and F = {/fc|l < k < np} is the list of triangles. We assume S 
is oriented and consistent. That is, neighboring triangles have the normals 
pointing to the same side of the surface. For a given vertex v in V, we say 
that another vertex w is a. neighbor vertex of v if there is a triangle f in F 
such that V and w are both in /. We denote by Ng{v) the set of neighbor 
vertices of v in y and m the number of points in Ng{v). We also denote 
by T{v) the set of triangles f in F with v E f. If the triangle / is in T{v), 
we say that / is incident to v. The area of a triangle / is denoted by |/|. 



3.1 The least square method of Chen and Schmitt 

Chen and Schmitt [lj(1992) provided an algorithm to estimate the principal 
curvatures by the Euler formula. Their main idea is to choose a suitable 
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coordinate system {ri, on the tangent space. Given a unit vector t in the 
tangent plane, the Euler formula (jHI) will give 

Kit) = cos\e + Oq) + sin2(0 + ^o) (7) 

where is the angle between the principal direction ei and ri and 9 is now 
the angle between the unit vector t and ri. This equation ((7j) can be rewritten 
as 



K{t) = Ci cos^ 9 + C2 cos^ e + Cs cos^ 9 (8) 

for some constant Ci, C2 and C3. 

Given a vertex Vi in Ng{v), we can obtain an unit vector tj in the tangent 
plane of S at the vertex v by 

^ (fj-t;)- <Vi-v,N > N 
' \\{v,-v)-<v,-v,N>N\\ 

where N is the normal vector at v. The normal curvature kn{t) along the 
unit vector ti can now be approximated by 

, , , 2 < Vj - V, N > , . 

KiU) = ^ 'rr. . 10 

\\v — 

See Chen and Wu |2| (2004) for the discussions. Chen and Schmitt used the 
least square method to find the constants Ci, C2 and C3: 

min I (Ci cos^ 9i + C2 cos9i sin 9i + C3 sin^ 9i) — A:„(tj) |^ (11) 

i 

where 9i is the angle between ti and ri. The principal curvatures can 
then be solved from the constants Ci, C2 and C3 via the following relations: 



K,i cos^ ^0 + '^2 sin^ ^0 = 

2(ki — K2) cos6'o sin6'o = C2 (12) 
K2 cos^ 9o + Ki sin^ 6*0 = C3 

Moreover, the principal directions can also be computed from 



= cos(-6'o)ri + sin(-6'o)r2 
V2 = sin(-6'o)ri + cos(-6'o)r2 
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3.2 The integral method of Taubin 

To find the principal curvatures ki and ^2, Taubin [TU](1995) considered the 
following integral of a symmetric 3x3 matrix: 

^ = IJ Uv>v^d9 (14) 

where v = vi cos 6 + V2 sin 9. Taubin showed that the matrix B can be 
decomposed into 

B = E^(7 (15) 
V may/ 

where E = (^Vi V2) is the 3x2 matrix given by the principal directions. 
Moreover, the principal curvatures are given by the relations: 

( Ki = 3mi - 1712 

\ K2 = 3m2 — mi 
Taubin approximated the integral (jl4j) by the finite sum 

m 

B = Y,^ikn{ti)titf, (17) 

i=l 

where the unit tangent vector tj is given in © and the normal curvature 
is approximated as in equation (jTU)). The weight Ui is chosen to be propor- 
tional to the sum of the areas of the triangles incident to both v and Vi with 
Yl^i'^i — 1- Taubin decomposed the matrix B with a suitable transforma- 
tion and a rotation. In (Chen and Wu [2]), the centroid weights were used 
for Equation (fTTj) and they gave more accurate results. 

3.3 Related works 

Flynn and Jain [5^(1989) used a suitable sphere passing through four vertices 
to estimate curvatures. Meek and Walton j7j(2000) examined several meth- 
ods and compared them with the discretization and interpolation method. 
Gatzke and Grim[n](2003) systematically analyzed the results of computation 
of curvatures of surfaces represented by triangular meshes and recommended 
the surface fitting methods. See also Petitjean |9j(2002) for the surface fitting 



6 



methods. Meyer et al. [S] (2003) employed the Gauss-Bonnet theorem to esti- 
mate the Gaussian curvatures and introduced the Laplace-Beltrami operator 
to approximate the mean curvature. 

4 The differential approach: our method 

In this section we shall describe a new, simple and effective method to approx- 
imate the differential of the Gauss map. In order to simplify the presentation, 
we shall follow the ideas developed in Wu, Chen and Chi|Bj(2005) about the 
gradient and the Laplacian of a function defined on a triangular mesh. 

Consider a triangular mesh S = (y,F), where V = {vi\l < i < n^} is 
the list of vertices and F = {/fc|l < i < np} is the list of triangles. Let 
g he a. function on V. First, we can extend the function g to a piecewise 
linear function, still denoted by g, on S as follows. Given a face f in F with 
vertices Vi, Vj and Vk, every point p in f can be written as an unique linear 
combination of Vi, Vj and Vk- That is, p = avi + (3vj + jVk with a, 7 > 
and a + (3 + J = 1. Then we define g{p) by 

g{p) = ag{vi) + (3g{v.j) + 7^(ffc) (18) 

Thus, the function g is affine on each face f of S and hence is differentiable 
on /. The gradient (Vg) f of g on the face / at the vertex Vi can be computed 
from 

C^g) f{vi) = a{vj - Vi) + h{vk - Vi) (19) 
where the coefficients a and h are determined by the relations 

g{vj) - g{vi) =< {Vg)f{vi),Vj - v, >, .^q) 
9{vk) - givi) =< (yg)f{vi),Vk -Vi>, 

A direct computation gives 

f(^\ ^ [< - Vi,Vj - Vi >, < Vj - Vi, Vk-Vi>\ V g{vj) - g{vi)\ 
\bj \<Vj -Vi,Vk-Vi>, <Vk-Vi,Vk-Vi>J \g{vk) - g{vi) J 

(21) 

To obtain the gradient V g{y.i) of g on S at the vertex Vj, we use the 
weighted combination method. Namely, we set 

Vg{v,)= J2 ^f[i^9)M)] (22) 

/6T(i,,) 
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with u}f>0 and E/gt^) ^/ = 1- According to Chen and Wu 001(2004, 
2005), we shall use the centroid weight for the gradient formula ()22|) . The 
centroid weights are 



-/ = ^ (23) 

^/eTK) \\Gf-v,\\^ 

where the centroid Gk of the triangle face / is determined by 

= !^l±^!l±^ (24) 
3 

Next we consider the differential dN of the Gauss map A^. As in Chen 
and Wu[2](2004), we define the normal vector N{v) at each vertex in F by 
the centroid weights: 

Niv) = ^f^n^^^f^^ (25) 



where Nf is the unit normal to the triangle face / and the weight cu/ is given 
in (j23|l . The normal vector N{v) has three components: 

N{v) = {n,{v),n2{v),ns{v)f (26) 

and its components ni{v) are functions on V. Thus, we can compute their 
gradients and obtain the differential dNy of at the vertex v as 

dN, = (Vni(v), Vn2{v),Vn3{v)) . (27) 

Note that the differential dN, is a 3 x 3 matrix. In order to obtain 
the linear map dN, from the tangent space TS, into itself, we choose an 
orthonormal basis {ci, 62} for the tangent space TSp. That is, the vectors ei, 
62 and N{v) form an 3-dimensional orthonormal basis. Then the differential 
dNy can be realized by a 2 x 2 matrix, still denoted by dN, = {ciij) , and the 
entry Ojj is given by 

ttij =< ei,dN,{ej) > . (28) 

Therefore, we can estimate the Gaussian curvature K and the mean curvature 
Hhy 

K = aiia22 — 012021 /r,q\ 

TT _ aii+a22 
2 

The principal curvatures and the principal directions can also be com- 
puted from the eigenvalues and eigenvectors of dN,. 
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5 Computational Results 



Taubin's method to estimate the tensors product of curvature on triangular 
mesh is very useful in CAD. Furthermore, Chen and Wu jT] (2004) provided 
a better choice in Taubin's algorithm by centroid weights. We will compare 
Taubin's method by centroid weights and our method. In our tests, we 
consider the random polynomial surface, 

i: = {{u,vj{u,v))} (30) 

with 



CijU V 



i=Q j=0 

We compute the Gaussian curvature at the vertex p = (0, 0, /(O, 0)) on some 
different random surfaces. The set of neighbors of p is constructed by 

{{riCos9i,riSm9i, f{ri cos 9i,ri sin 9i))\i G {1,2,-- - ,nv}} (31) 

In equations and EH C-ij IS 8b random number in the interval [—5,5], 
{9o, 9i, - ■ ■ 9n^} is a random partition of [0, 2n] such that < \9j^i—9j\ < I.Qtt 
for each j and m,n,ri,nv G M"*" are some random positive values. And we 
estimate the error of Gaussian curvatures by the formula 



Err{K) 



\K-K,\ 
K 



where K is the real Gaussian curvature at vertex v and is the Gaussian 
curvature at vertex v obtained by Taubin's method or our method. 

From figure HI to El we tests 1,000 random surfaces. For each surface, we 
compute the error of average of 10,000 different random partitions. From 
these figures, our method is better than the Taubin's method. In Figures 
m El and El we test the effect of different partitions. For each partition, we 
choose 10,000 differential random surfaces and estimate the error of average 
and standard derivation. Obviously, Our method is more stable than the 
Taubin's method. 

Final Remark: 

This method is conceptually simple and more natural than the normal curva- 
ture method of Chen and Schmitt and Taubin. In the next section, we shall 
show that this method also yields more accurate results. In Wu, Chen and 
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Figure 1: The effect of different surface of Taubin's metliod 



Chi|llj(2005). tlie autliors also develop a differential theory for triangular 
meshes. The gradient, Laplace-Beltrami operators are discussed. Moreover, 
this method also works for boundary vertices and for polyhedron meshes. 
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